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Introduction 

The scope of the work in this presentation focuses on the development of algorithms for 
the integration of rate dependent constitutive equations. In view of their robustness; i.e., 
their superior stability and convergence properties for isotropic and anisotropic coupled 
viscoplastic-damage models, implicit integration schemes have been selected. This is 
the simplest in its class and is one of the most widely used implicit integrators at present. 

Viscoplastic Models 

Several viscoplastic models have been proposed and developed to treat the complex 
time dependent viscoplastic behavior of metals, alloys and composites at high tempera- 
ture. The deformation behavior of materials at high temperature involves energy dissipa- 
tion and material stiffness variations due to physical changes in the material's 
microstructure. Consequently, thermodynamic arguments have often been utilized as a 
foundation on which phenomenological constitutive laws may be formulated. The com- 
plete potential-based class of inelastic constitutive models exhibit a number of unique 
advantages from both a theoretical and a computational standpoint, for example, the 
symmetry of the resulting consistent tangent stiffness matrix, and possesses a form 
which is convenient for further development of new deformation and damage models. 
The Generalized Viscoplastic with Potential Structure (GVIPS) [1] model possesses 
both the thermodynamic potential (Gibb’s function) and the dissipation function ( Q 
form). 

Another class of constitutive models are the Non-Assodative Viscoplastic (NAV) mod- 
els. the NAV models refer to those that have a partially (e.g. Q form only) or totally 
incomplete potential form. An example of a NAV model is that of Freed [2]. Recent work 
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has demonstrated that the above models may be modified to restore the complete poten- 
tial structure. 

For these two different classes of models, a general computational framework suitable 
for implementation of both is needed. 

Integration Schemes 

Computational algorithms for the integration of constitutive relations play a key role in the 
inelastic finite element analysis of engineering structures. Consequently, much research 
effort has been devoted over the years to the development and critical assessment of 
integration schemes for the rate equations in these material models. 

In early applications the explicit integration schemes, (i.e., forward Euler method) were 
predominate because of their ease of implementation, and because they do not require 
evaluating and inverting a Jacobian matrix. However, explicit integrators may not be effi- 
cient. That is, too many iteration steps may be required and convergence (stability) can- 
not be guaranteed. As a result, several alternative approaches have been used, for 
example, Gear's multi-step method [3] and Walker’s asymptotic method [4]. Note that 
every integration scheme has its own particular application domain and is problem 
dependent. 

The majority of recent work has emphasized the use of implicit integration methods in 
view of their stability and convergence properties. Based on the fully implicit, backward 
Euler scheme, the corresponding algorithmic (consistent) tangent stiffness arrays are 
derived from the integration rule, which are important for finite element solutions using 
(global) Newton-Raphson iterative methods. 

Line Search 

Although the implicit scheme is unconditionally stable, its successful application may still 
require proper selection of the size of the steps utilized. In this regard two factors are 
important: (i) accuracy, and (ii) convergence of the local iterations. A simple time subin- 
crementing strategy was found to be effective in obtaining accurate results especially 
when dealing with regions of discontinuity in the state space. However this was found to 
be insufficient to obtain a computationally efficient solution for a highly nonlinear problem 
such as viscoplasticity. When a large time-step size is chosen, too many subincrement- 
ing are needed, which leads to inefficiency. Thus a more sophisticated solution proce- 
dure, namely, a line search algorithm, is required to produce an effective robust solution 
algorithm. 

It is well know that classical Newton-Raphson is fast and stable only when the trial solu- 
tion is close to the converged value. Thus, the purpose of the line search algorithm is to 
guide the solution towards convergence by searching for a scalar multiplier that adjusts 
the amount of the increment vector to be updated within each iteration [5]. The concept 
of line search may be applied at either the global (structural) iteration level or at the local 
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(constitutive) iteration level. At the global level, the concept of the line search algorithm 
pertains to minimizing the total potential energy, that is, the work done by the residual 
force due to the iterative displacement. It has been suggested that the line search be 
incorporated with a consistent tangent stiffness and that the use of the line search is 
essential for robust performance of Newton’s method [6,7]. It also demonstrated that in 
elasto-plastic analysis convergence is not guaranteed unless the global line search is 
used [8]. On the local (constitutive) level, line search is used to adjust the suitable incre- 
ment of stress and internal variables to guarantee the convergence of the local iterations. 
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MOTIVATION 


Two major obstacles for fully utilizing recent time-dependent/hereditary 
constitutive models in practical engineering analysis 

• Lack of efficient and robust integration algorithms 

- Coupled system of tensorial rate (differential) equations or general kernal 
convolution integration 

- Increased mathematical complexity and associated numerical stiffness 

• Difficulties associated with characterizing large number of required material 
parameters 


Ro-i 


OBJECTIVE 


• Develop a robust and efficient integration algorithm for viscoplastic 
constitutive equations 


Fig. 2 
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BACKGROUND 

• Computational algorithms for integrating the constitutive models are a key 
component for an efficient inelastic finite element analysis 

• Two classes of integration schemes are: iterative and non-iterative 

• Iterative 

- Fully implicit scheme 

- Requires local iterations 

• Non-iterative 

- Semi-implicit 

- Fully explicit 

- Nq local iterations 

- Usually less “overhead" 

• History dependent integral representation 

- Full history data storage 

Fig. 3 


INTEGRATORS 


General differential form: e n + j = e n + At 


[(l-a)e'+<*4 +1 ] 


• Fully explicit: (a = 0) Forward Euler 

• Fully implicit: (a = 1) Backward Euler, (a = 1/2) Midpoint rule 

- e must be evaluated at n + a, which requires a local iterative procedure 

• Semi-implicit: (0 £ a £ 1 ) 

- to avoid iterations, some methods employ an approximation for e n + j 


- these are referred to as toward gradient methods, i.e., £ n+ i is approxi- 
mated in terms of known quantities at time n using a Taylor series expan- 
sion, 


./ 

e n + 1 


/ 

= E n + 


da 


Ao n + 


n 




Fig. 4 
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INTEGRATORS 


General Integral Form: 

t 

A = \K(t, T)e(T;)<rc 
0 

A = Typical state variable 
K(t, x) = Kernal function, e.g., 

• exponential 

• power form 

• functional derivative 


Fig. 5 


MODEL CLASSES INVESTIGATED 

• Differential/Internal variable type 

- nonassociative/dynamic recovery (NAV) 

- fully associative/nonlinear kinematic hardening (GVIPS) 


Fig. 6 
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NON-ASSOCIATIVE MODELS (NAV) 

FREED/WALKER VISCOPLASTIC MODEL (1993) 
(DYNAMIC RECOVERY) 

Flow Law: e = f(J,D)T T = M(g-g) a = a s + a l 


fU ' D) = M irj smh ^" 

M = isotropic/deviatoric tensor operator 
Evolution Laws: 

a s = 2 Z[tf / - « A ] a, = 2z[«/ - *,*,] 

5 / = ¥«/ D = qj-q D 

a = “short-term” back stress; a, = “long-term” back stress 


Fig. 7 


POTENTIAL BASED MODELS 

GENERALIZED VISCOPLASTICITY WITH POTENTIAL STRUCTURE, GVIPS 
(NONLINEAR KINEMATIC HARDENING) 


Gibb’s Potential: »F(g, g) = ¥ e (g) + 'P'(g) Dissipation Potential: Cl = G(g. g) 
Flow Law: z = ^ * /(F)r r = M(g-g) 

M = isotropic/anisotropic tensor operator 


/(F) = 


2p 


,-l8Q d'PdQ 
Evolution Law: g = L = 


= 4 _1 ^e 7 - n = Mg 


L _1 = 


da dadada 

= hardening stiffness operator 


r d 2 'P‘ l 1 


3a9a 


Fig. 8 
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IMPLICIT INTEGRATOR 


~n + 1 ~ ^n + ^-n + 1 
• £ , i is the vector of state variables 



f<V 



For the NAV model: Z = 

^sn 

and for GVIPS: Z = 

fV 

~n 

*ln 

~ n 



y d J 




o stress 
~ n 


g n back stress (s and l denote short and long term) 
D drag stress 


Fig. 9 


IMPLICIT INTEGRATOR 

dZ n + 1 is the increment in state variables and may be expressed as: 

«» + l - (*£>'' l *« + l 

- where is the iterative Jacobian matrix of state variables 

'- n +l is the residual function of state variables, e.g., # n + l = 

- and the residual functions for NAV are 

^ = 5„ + r?„T e (Ae + A//r fl + 1 ) 


V 


~a, 


V r d) 


R„ = a , 
a s,l ~ s ’ 1 


" «c / " 2AtH S ifZZn + 1 + 2Atg S l~~S l 

n+l s,l n S ’ 1 ” + 1 S ’ 1 s,l n + 1 


R = R +R 
~a ~a -a, 
s l 


R D ~ D n + 1 " D n ~ At ^ q J ~ q D* 
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LINE SEARCH 


• The factor ti is a scalar (0 s ti £ 1 ) that adjusts the step size to optimize the 
iterative solution 

- ti is obtained by a ling search algorithm 

I k k 

- the objective is to minimize the dot product s ( = U?„ + j + i 


ratio s i 
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VALIDATION TESTS 

CYCLIC CREEP RELAXATION 
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Inelastic Strain 


CYCLIC TEST RESULTS 

NAV MODEL 



CREEP AND RELAXATION RESULTS 

GVIPS MODEL 

CREEP RELAXATION 



Fig. 14 
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CREEP AND RELAXATION RESULTS 

NAV MODEL 

CREEP RELAXATION 
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NONPROPORTIONAL LOAD 

NAV MODEL 

Nonproportional loading path Nonproportional loading 

results for NAV (121 °C) 



Fig. 16 
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LINE SEARCH COMPARISONS 

CYCLIC TEST (NAV MODEL) 
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CREEP TEST (GVIPS MODEL) 
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RELAXATION TEST (GVIPS MODEL) 
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SUMMARY/CONCLUSIONS 

• Implicit integration algorithm provides unconditional stability 

- for both the creep and relaxation tests the explicit integrator (forward Euler) 
failed at 1 00 steps 

- implicit succeed using only 2 steps 

- “large” time increments for an efficient solution (computation time savings) 

• Accuracy is consistent with first order formulation 

- creep: 2 steps, 8% error relaxation: 2 steps, 2% error 

• Important for life prediction studies that require many analysis load cycles and 
an efficient integrator 

• Current algorithm used in material parameter estimator 

- analysis is performed repeatedly during optimization: requires efficiency 

- parameters vary during optimization: requires robust integrator 

Ftg. 18 
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FUTURE RESEARCH 

• Currently, implicit integration algorithm has been successfully implemented into 
the MARC user subroutine HYPELA 

• Next: Implement algorithm into ABAQUS user subroutine UMAT 

• Organize computer code to allow easier implementation of new constitutive 
models 

• Develop for combined differential/integral hereditary representations 


Fig. 19 
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